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Percus showed that approximate theories for the structure of nonuniform hard sphere fluids can 
be generated by linear truncations of functional expansions of the nonuniform density p(r) about 
that of an appropriately chosen uniform system. We consider the most general such truncation, 
which we refer to as the shifted linear response (SLR) equation, where the density response p(r) to 
an external field <f>(r) is expanded to linear order at each r about a different uniform system with a 
locally shifted chemical potential. Special cases include the Percus- Yevick (PY) approximation for 
nonuniform fluids, with no shift of the chemical potential, and the hydrostatic linear response (HLR) 
equation, where the chemical potential is shifted by the local value of cj>(r). The HLR equation gives 
exact results for very slowly varying cf>(r) and reduces to the PY approximation for hard core 0(r), 
where generally accurate results are found. We show that a truncated expansion about the bulk 
density (the PY approximation) also gives exact results for localized fields that are nonzero only in 
a "tiny" region whose volume V* can accommodate at most one particle. The SLR equation can 
also exactly describe a limit where the fluid is confined by hard walls to a very narrow slit. This 
limit can be related to the localized field limit by a simple shift of the chemical potential, leading to 
an expansion about the ideal gas. We then try to develop a systematic way of choosing an optimal 
local shift in the SLR equation for general <f>(r) by requiring that the predicted p(r) is insensitive 
to small variations about the appropriate local shift, a property that the exact expansion to all 
orders would obey. The resulting insensitivity criterion (IC) gives a theory that reduces to the HLR 
equation for slowly varying 4>(v), and is much more accurate than HLR both for very narrow slits, 
where the IC agrees with exact results, and for fields confined to "tiny" regions, where the IC gives 
very accurate (but not exact) results. However, the IC is significantly less accurate than the PY and 
HLR equations for single hard core fields. Only a small change in the predicted reference density is 
needed to correct this remaining limit. 



I. INTRODUCTION 

In this paper we discuss approximate methods for de- 
termining the structure and thermodynamics of nonuni- 
form hard sphere fluids in a general external field 4>(r). 
The field can directly describe the effects of fixed solute 
particles, confining walls, and other sources of nonunifor- 
mity on the hard sphere flui d- 1 '"'S' and the model serves as 
a reference system through which properties of nonuni- 
form simple liquids can be determined by density func- 
tional theory 4 or molecular field theory^ 

We start from a general theoretical perspective first de- 
scribed in a classic article by Percus^ and reviewed below 
in section [H] Percus suggested that approximate theo- 
ries for the density response to (j>(r) could be generated 
by truncations of formally exact functional Taylor series 
expansions of the nonuniform density p(r) about that of 
an appropriately chosen uniform system, fn particular 
he argued that the density response to a hard core <j>(r) 
could be accurately described by a simple linear trunca- 
tion of an expansion about the bulk density. This yields 
the Percus- Yevick (PY) approximation for nonuniform 
fluids,' discussed in section ITTTl and indeed this gives 
rather accurate results for hard core solutes of varying 
diametersi&iLiSiAi The linear truncation is the key to de- 
veloping a practically useful theory, since higher order 



terms in the expansions are very complicated. Unfortu- 
nately, the results of the PY approximation deteriorate 
very quickly for more extended fields, especially ones with 
attractive regions, where the linear extrapolation from 
the bulk becomes very inaccurate. 

In sections IIVI and [V] we exploit an exact property of 
the grand canonical ensemble, where both the chemical 
potential fi B and external field 4>(r) can be shifted by a 
constant without changing any properties^ This allows 
us to consider a more flexible expansion where the density 
response p(r) at each r is expanded to linear order about 
a different uniform system with a locally shifted chem- 
ical potentialAifliii We can use the additional freedom 
permitted by this local shift to minimize errors arising 
from the linear truncation of the expansion. Indeed, we 
show below that there are several different limits where 
exact results can be obtained from a linearly truncated 
expansion by using an appropriate shift. 

When the field <f>(r) is very slowly varying, it is natu- 
ral to expand about a uniform system where the chem- 
ical potential is shifted by the local value of the field. 
This leads to the hydrostatic linear response (HLR) equa- 
tion, which we derived earlier using physically motived 



arguments. 



The HLR equation gives exact results 



for very slowly varying </>(r) and turns out also to give the 
same predictions as the PY approximation for hard core 
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</>(r). Thus it can be successfully used for a much wider 
class of fields (both repulsive and attractive). However, 
the specific choice for the shift made in the HLR equation 
is based on the local value of the field, and this can give 
poor results in certain special cases, particularly when 
there are rapid variations in <f>(r) in a confined region 
of space much smaller than the correlation length of the 
bulk fluid. 

But we show in section that in two such cases a lin- 
early truncated expansion can again give exact results 
with different choices of the locally shifted system. In 
particular we show that exact results can be found by 
expanding about the bulk density (i.e., the PY approxi- 
mation) for external fields 4>{v) that are nonzero only in 
a tiny region of space whose volume is so small that 
the center of at most one hard sphere can be accommo- 
dated within it. The fact that the PY approximation can 
give exact results for such confined or "tiny fields" does 
not seem to be known in the literature. This limit is of 
more than academic interest since some recent advances 
in density functional theory have resulted from requiring 
that it be satisfied^ Moreover, we show how expansions 
used in the tiny field theory can be easily modified to ex- 
actly describe a seemingly very different problem, a nar- 
row slit, where the potential is infinite except in a narrow 
rectangular slit of width L s . In the limit as L s — > 0, we 
show that a truncated expansion about the ideal gas can 
give exact results^ 3 - whereas the PY (and HLR) equation 
gives very poor results in this limits 

The fact that different shifted linear truncations can 
exactly satisfy all these distinct limits suggests that a 
generally useful theory might arise from finding a local 
shift that in some sense minimizes the errors arising from 
the linear truncation. In section IVT1 we present our first 
attempt along these lines, where the local shift is chosen 
self-consistently so that the predicted density response 
p(r) is insensitive to small variations about the proper 
local choice, a property that an exact expansion to all 
orders would satisfy. This "insensitivity criterion" (IC) 
effectively generates an expansion about a local uniform 
system whose density can be interpreted as a smoothed 
weighted average of the full density p(r), reminiscent of 
results in certain versions of weighted density functional 
theoryAi^ 

Results of the various theories for the density response 
to several different fields are presented in section IVIII 
In general the IC performs very well: it reduces to the 
exact HLR equation for slowly varying fields, gives exact 
results in the narrow slit limit, and generally accurate 
(though not exact) results for tiny fields. However, the 
IC approach is significantly less accurate than the PY 
or HLR equations for a single hard core field, where it 
substantially overestimates the contact density at high 
density. Since this is one limit where the PY and HLR 
equations do give satisfactory results, this may not be of 
much practical importance. But it would be better to 
have a general approach that remains accurate for these 
cases as well. Suggestions for further work along these 



lines and final remarks are given in section lYUH 



II. EXACT EXPANSIONS FOR THE DENSITY 
RESPONSE TO AN EXTERNAL FIELD 

The potential distribution theorem of Widom ;L — relates 
the density p(r) of a single component hard sphere fluid 
in an external field 0(r) to the probability P(r; [p]) of 
inserting a hard sphere test particle at the position r: 

p(r) = A- 3 e /3 ^ B -^«]p(r; [p]). (1) 

The direct influence of the external field on the test par- 
ticle is excluded in determining P(r; [p]), which depends 
only on the intermolecular interaction energy between 
the test particle and the hard sphere fluid. Equivalently, 
P(r; [p]) is the probability that a cavity whose radius is 
equal or greater than the diameter d of the hard sphere 
particles exists at the particular point r, since only then 
can the test particle be inserted. Here p B is the chemical 
potential, A the thermal de Broglie wavelength, and [p] 
specifies the functional dependence on the density pro- 
file. P(r; [p]) can be formally reexpressed^ 3 - in terms of 
the one-body direct correlation function c^\r; [p]) as 

P(r;[p]) = e c<1, ( r ^. (2) 

By expanding P(r±; [p]) at a given ri in a functional 
Taylor series about a uniform fluid at some density p, we 
arrive at a formal expansion for the nonuniform density 
directly related to an expansion suggested by Percus&i& 

p(n) = A - 3 e^ B -^ ri )]+ c(1) ^ 

x{l + J dr 2 c^(r 12 -p)(p(v 2 )-p) 
+Q(r x ) + ...}, (3) 
where the quadratic term in the expansion is 

Q(ri) = J dr 2 J dr 3 {cS 2 \r 12 ; p)cS 2 \r 13 ; p) 

W 3) (ri,r 2 ,r 3 ;p)} (p(r 2 ) - p)(p(r 3 ) - p)(4) 

The c^(ri,r2, r„; [p]) are defined by successive func- 
tional derivatives of cW(r; [p]) with respect to the singlet 
density, i.e., 

c(")(, 1 ,r 2 ,...,r m;M )^ fc(n " 1) ^; , "; ,r "" i:M) - ( 5 ) 

op{r n ) 

Since we expand about a uniform fluid state, 
c( 2 )(ri,r 2 ; [p]) = c( 2) (r 12 ; p), where r 12 = |r 2 - ri|, due 
to translational invariance. Although the density p(r) 
can have discontinuities caused by the discontinuities in 
the external field <j>(r), representing, e.g., a hard wall 
or a spherical cavity, P(r; [p]) is always continuous and 
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smooth^ The expansion in eq |3] is designed to take ad- 
vantage of the smoothness of P(r; [p]). The hope is that 
with proper choice of p the expansion can be truncated 
at some low order and a relatively simple theory for p(r) 
will result. 

Percus considered many other expansions as well, some 
of which might seem even more promising & For example, 
by directly expanding c^\ti; [p]) in eq|2in a Taylor se- 
ries, we are guaranteed that the resulting approximation 
for the density after a truncation is always nonnegative, 
an exact and desirable property not always produced by 
truncations of eq [3] However as discussed in sections IIIII 
and Ivl there is a well defined limit where the expansion 
in eq 01 truncates exactly, in contrast to the correspond- 
ing expansion for c*- 1 ^ . Moreover we will show that both 
the HLR and the PY equations can be derived from eq 
|3| This suggests that it offers a versatile starting point 
for further research. 



is known to give reasonably accurate results for the den- 
sity response to larger spherical cavities (including the 
hard wall limit), where <f>(r) is infinitely repulsive inside 
a spherical region of radius R, but zero elsewhere. For 
such potentials the exact result p(r) = is trivially ob- 
tained from eq[S]in regions where <p(r) is infinite because 
of the Boltzmann factor e~^^ r \ 

Despite this success, the linear extrapolation using the 
bulk fluid in eq [5] would be expected to give poor re- 
sults for external fields that remain finite and vary over 
extended regions, especially in negative regions of the 
field where errors in the truncated series can be greatly 
magnified by the Boltzmann factor e~P^ r \ These are 
limits where the PY approximation is known to be very 
inaccurate ^ 



IV. EXACT SHIFTING PROPERTY OF THE 
GRAND ENSEMBLE 



III. CONFINED FIELDS AND THE PY 
APPROXIMATION 

One limit where the expansion in eq [31 is particularly 
useful is when the external field <fi is nonzero only a region 
much smaller than the correlation length of the fluid. 
The analyticity of P(r; [p]) then ensures that its values 
in the tiny region where <fi is nonzero can be accurately 
determined by making use of a low order extrapolation of 
its values outside, i.e., where <f> vanishes. For such locally 
confined fields, it seems clear that the expansion in eq|3| 
should be about p = p B = p{p B ), where p{p) gives the 
density of the uniform system with (f> = as a function 
of the chemical potential p. In terms of the quantities 
appearing in eqs^ and [21 this can also be written as 



P 



(6) 



It seems plausible that expansion only to linear order in 
eq[3]could then give an accurate description of the fluid's 
density response to a very confined field: 



A -3 e /3[/x B -#r 1 )J+c< 1 >(p B ) 

I + / dv 2 c^{r l2 ;p B ){p{v 2 )-p B ) 
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This qualitative argument will be made more precise in 
section Using eqEI this equation can be rewritten as 



1+ / dr 2C W(r 12 ;p B )(p(r 2 )~p B ) 



p( ri ) = p B e~^ 

(8) 

If the direct correlation function c* 2 >(ri 2 ; p B ) of the uni- 
form bulk fluid is known, cq[S]can be solved for the den- 
sity p(r) induced by the external field <j>(r). 

Equation [S] is the PY approximation for nonuniform 
fluids. 6 We argue below that it gives exact results for 
any sufficiently localized 4>(v), as suggested by the ex- 
trapolation argument. Moreover, the PY approximation 



When the linear truncation is inaccurate, it may seem 
difficult to make further progress, since the higher order 
terms in eq|3|are too complicated to use in practical cal- 
culations. However, as noted by Percus^ one does not 
have to expand eq|3 about the bulk density defined by eq 
[51 nor does one have to expand about the same bulk state 
for each r value of p(r). We can use this additional flex- 
ibility to greatly extend the accuracy of different linear 
truncations of eq[3] 

As is well known, it is the combination p B — <jf>(r) that 
determines the density profile in the grand canonical en- 
semble, and not p B and cf>(r) separately. 4 When both the 
chemical potential and the external field are shifted by 
the same constant, the system's properties should thus 
remain unchanged. This exact shifting property^ of the 
grand ensemble will play a key role in what follows. 

Using this shifting property in eq[21 we consider a gen- 
eral shifted chemical potential 



and a shifted external field 



(9) 



ri (r) = 0(r) - a r \ (10) 

both of which are shifted from the original p B and cf>(r) 
by a constant a ri that in principle can depend paramet- 
rically on the point ri about which the expansion is car- 
ried out, as indicated by the superscript ri. The shifted 
chemical potential p ri defines at each ri an associated 
bulk system with a uniform density 



P ri ^p(M ri ), 



(11) 



whose correlation functions are used in the expansion. 
Equation [3 thus becomes 

p(r x ) = A~ 3 e /3[ ' iri ^ ri(ri)1+c<1)( ' iri) 

x{l + J dr 2 c^(r 12 ;~p^ )(p(v 2 )-p^ ) 

+Q(ri) + ...} (12) 
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In principle (assuming convergence of the infinite se- 
ries), if we could accurately evaluate all terms in the Tay- 
lor series, the same exact result for p(i"i) would be found 
for any choice of a 1 " 1 due to the shifting property. But 
this is hopelessly complicated in general and approximate 
values for p(ri) from truncations of the series do depend 
on the choice of the a ri . 



V. SHIFTED LINEAR TRUNCATIONS 

Our strategy is to try to choose the a ri or the p ri 
at each ri in such a way that a self consistent solution 
for p(ji) based on a simple low order truncation of the 
series in eq|3]can give accurate results. In particular, we 
suppose that the a ri can be chosen by some argument to 
be specified later so that the expansion can be truncated 
at linear order to a good approximation. We thus arrive 
at a very general starting point, which we refer to as the 
shifted linear response (SLR) equation: 



p(r x ) = A~ 3 e' 3[ ' iri ^ ri(ri)1+C<1)( ' 5ri) 

1+ / dv 2 c^\r 12 -(f 1 ){p(v2)-p ri ) 



1+ / dr 2C ( 2 )(r 12 ;p ri )(p(r 2 )-p ri ) 



(13) 



Specific choices of the a Vl or p ri will lead to different 
approximations for p(i"i). The PY approximation arises 



from the choice a ri = or p r 



As we argued 



above this choice should give very accurate results for 
sufficiently localized fields and it is known to give a good 
description of the response to hard core solutes. 

However, when the external field is extended but slowly 
varying, a different choice of a Tl is clearly more appro- 
priate. As discussed in detail in References Wi fiol and ITU 
the HLR equation is very accurate in such cases. This 
uses the hydrostatic shift, where the external field is lo- 
cally shifted at each ri to be zero. This corresponds to 
the choice 



= 0(ri 



(14) 



in eqED S0 that (j> Tl (ri) = 0. The associated uniform 
density given by the shifted chemical potential for this 
particular choice of a ri is denoted p ri and it satisfies 



p^ =p(p B -^(r 1 )). 



(15) 



Unlike the PY approximation, the HLR equation 
builds in the invariance of the grand canonical ensem- 
ble when both the chemical potential and the external 
field are shifted by a constant. It is an excellent approx- 
imation when the external field is slowly varying, and 
the expansion in ea 1121 about the hydrostatic density eq 
^Jconverges rapidly^ Moreover the HLR equation gives 



the same results as the PY equation for hard core fields. 
In general the HLR approximation is much more accurate 
than the PY approximation for a wide range of external 
fieldsiMi 

The SLR ea uat ion 1131 with the particular choice a ri = 
</>(i"i) thus provides an alternative derivation of the HLR 
equation, in addition to the previous physical arguments 
based on optimizing linear response by use of the hydro- 
static shifti V ! 11 However, as discussed below, there are 
some limits where the HLR choice for a Tl gives poor re- 
sults. The SLR equation permits other choices for a ri 
and thus provides additional flexibility that can lead to 
improvements of the HLR equation. With proper choice 
of a r , eq^Jcan be used to describe and bridge several 
limits that both the PY and the HLR equations fail to 
capture. These limiting cases will be addressed in the 
following. 



A. Tiny Fields 

Consider first a very localized field 0(r) that is non- 
vanishing only within a spherical volume of radius 
d/2, with d the hard sphere diameter: 



(r) + 0, |r| < - 
= 0, otherwise. 



(16) 



This volume is so small that it can simultaneously accom- 
modate the centers of at most one hard sphere particle. 
We call localized fields that are nonzero only in such a 
tiny region tiny fields. A special case is a tiny cavity, 
where <p is infinite in . The density response of a hard 
sphere fluid to any tiny field can be determined exactly, 
as we now show. 

We start with the grand partition function S[</>] when 
a general external field (j> is present in the fluid: 



oo N 
Z 



N=0 



-0 4>(n)-l3V N ({r}) 



oo N 
Z 



N=0 



Here Z/v[0] is the canonical partition in the presence of 
the field, ^at[0] refers to that of the uniform fluid when 
the external field is zero, z = exp(f3p B )/A 3 , and V/v({r}) 
is the intermolecular interaction potential between the N 
fluid particles. Introducing the Mayer /-function for the 
external field 



/V) 



= „-/3<«r) 



1. 



(18) 



eq 1171 can be rewritten in terms of an expansion about 
the uniform fluid with d> = 0: 
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00 Z N N r 



(JV-n)!n! 
/ dr 1 dr a ...dr„/*(r 1 )/*(r 4 ).../*(r n ) £ 



/*(ri)/*(r 2 ).../*(r„) / dr„ +1 ,...,dr w e-' 3y "« r » 
dr n+1 ,...,dr Are -' 3y "« r » 



JV>n 



(iV-n)! 



oo - 

S ME- / dr 1 dr 2 ...dr n /*(r 1 )/^(r 2 ).../^(r n )p(»)(r 1 ,...,r n ;[0]) 



(19) 



n=0 



Here p( n )(ri, r„; [0]) is the n-particle distribution func- 
tion in the uniform grand canonical ensemble with chem- 
ical potential p B . 

For general extended fields this formal expansion does 
not converge rapidly. However it can be very useful when 
the field is confined to a small localized region of space 
since /^(r) is nonzero only where 4> is nonzero. In partic- 
ular for tiny fields the expansion ea 1 191 must truncate ex- 
actly due to the vanishing of the p^ n '(r±, r n ; [0]) when 
more than one hard particle is simultaneously within the 
volume V^. By functional differentiation of eq^Dwe can 
also find an expansion for p(r; [</>]) that similarly trun- 
cates. Thus we have exactly for tiny fields 



1 + / deip»f* (n) 



(20) 



and 

p(n;M) 



1 + / dv 2 p B f*{v 2 ) 

I - / dv 2 p B g^{r 12 ;p B )^{v 2 ) 



(21) 



The exact result for the density response to a tiny cav- 
ity was originally derived in a different way by Reiss and 
CasbergpiS, Here g^(ri 2 ; p B ) is the exact radial distribu- 
tion for the uniform hard sphere fluid. This same formula 
clearly holds for a more general model with longer ranged 
pair interactions outside the hard core if the appropriate 
is used. 

As would be expected by the appearance of p B in this 
equation, one can show that the PY approximation eq|S| 
is consistent with this exact result for any tiny field. See 
reference for a straightforward but tedious derivation. 
In accord with the qualitative argument above, the linear 
extrapolation from the bulk into the tiny region is exact 
in this case. This corresponds in the SLR equation to 
the choice a r = 0, i.e., p r = p B and shows that the series 
indeed truncates exactly in this special case. 

However, a closely related limit highlights a general 
problem with the PY approximation. Consider a field 
4> c (r) that is a non-zero constant c outside a tiny region 
and any value 4>(r) inside. This can immediately be 
shifted to be the type for which eqs |20 an d |^ hold by 



r 



making the choice a r — c. For such a field (j) c we then 
have 



p(n;M) 



1 + / dr 2/ 5 c / c (r 2 ) 
x 1+ / dr 2 p c g^(r 12 ;p c )f c (r 2 ) 



(22) 



Here p c = p(p B — c) and f c is defined as in ea 1181 with 
the shifted tiny field 4> c (r) — <fi c {r) — c. 

For a perturbation that varies significantly only in a 
very local region (compared with the correlation length 
of the particles), a particle situated in the perturbed re- 
gion, though directly affected by the field, should screen 
the perturbation from the rest of the fluid. The fluid's re- 
sponse thus essentially remains that of the uniform fluid 
outside the local region. The shifted field represents such 
a localized perturbation. Thus choosing the density p c 
to be that of the bulk environment, i.e., p(p B — c), will 
truncate eqE3 a t l° w order, leading to eq[22 And again, 
the choice a r = c in the SLR equation ^| gives this exact 
result. 

In this essentially equivalent case however, the PY ap- 
proximation, which always uses the unshifted p B , will give 
an incorrect result, even though it can exactly describe 
the tiny field limit when c = 0. Unlike the SLR equa- 
tion, the PY approximation does not build in the exact 
shifting feature of the grand canonical ensemble. This 
can cause significant errors for extended slowly varying 
fields. 

However, the HLR approximation uses the hydrostatic 
shift eq El an d thus will correctly describe the shifted 
bulk density p c in this case. Moreover it is exact for tiny 
cavities (tiny fields that are infinitely repulsive inside V^) 
since any finite value for a r inside the cavity will still give 
the correct zero density. However the HLR equation is 
not exact for more general finite tiny fields. The HLR ref- 
erence density p r would follow the variations in <f> inside 
the tiny region, contrary to the exact result with a con- 
stant p c everywhere. Rapidly varying tiny fields can thus 
generate noticeable errors in the HLR approximation, as 
will be shown in a later section where computational re- 
sults are reported. 
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B. Narrow Slits 

Another application that may at first seem to be very 
different from the tiny field case is when an extended 
external field confines the fluid to a region of reduced 
dimensions. For example, consider a hard sphere fluid 
confined between two planar hard walls forming a nar- 
row slit. We can get exact results for this case from eq 
by exploiting the shifting property of the grand en- 
semble. The confining field can be taken to be the limit 
of piecewise constant potentials defined so that 

4> c (z) = 0, < z < L s 

— c, otherwise, (23) 

in the limit where c — > oo. In this limit the fluid's den- 
sity will be zero except in the narrow region between the 
walls. Here L s is the effective width of the slit as seen by 
the centers of the fluid particles. 

If we formally introduce the uniform shift a z = c we 
have 

(f> c (z) = -c, < z < L s 

= 0, otherwise (24) 

and the shifted external field <jf{z) is non-zero only in the 
narrow slit region, similar to a tiny field. However, the 
shifted slit field is not strictly a tiny field as defined in the 
previous section where the expansion exactly truncates, 
because even when L s — > 0, many particles in principle 
can still be found in the slit, aligned in a two-dimensional 
layer along the walls of the slit. But once this shift has 
been made, the expansion in eqs 1191 and 1221 converges 
rapidly for small L s since the contributions from the in- 
tegration over the f c tend to zero. The shifted chemical 
potential fl c = pP — c tends to — oo, corresponding to 
an expansion about the ideal gas limit where the shifted 
bulk density p c tends to zero and cW(/5 c ) = 0. 

Inside the slit where p c — </> c (i"i) = p B — </> c (ri) = p B , 
we have a finite limiting density as L s — > given by 

^-^(n) = A -3 e /3[A c -0 c (ri)]+ c (1) (p c ) = A" 3 e /3/iB 

Equation 1221 then gives the first two terms in an ex- 
act (but non-truncating) virial-like expansion valid for 
narrow slits. Higher order terms can be determined 
straightforwardly from eq 1191 In agreement with pre- 
vious work^ there is a constant limiting lowest order 
density profile in the narrow slit given by 

p(z) = A- 3 e^ B , 0<z<L s 

= 0, otherwise. (26) 

For a value of p B corresponding to a dense uniform hard 
sphere fluid, this yields a very large limiting value for 
the reduced density in the slit p 3D d 3 = A _3 e /3M d 3 3> 1. 
However the density of fluid particles per unit area of 
the wall p 2D d 2 = A _3 e /3M d 2 L s tends to zero as L s — > 



0. Thus particles in the narrow slit are very far apart 
laterally and an expansion about the ideal gas limit is 
physically appropriate. 

Clearly the SLR equation can reproduce these exact 
limiting results if the proper choice a z = c (where c — ► 
oo) for all z is made. However, the PY approximation 
uses a z = everywhere in the SLR equation, while the 
HLR equation assumes a z — c outside the slit but a z = 
inside, and hence they both give incorrect results in 
the limit L s — > 0. Both theories correctly predict zero 
density outside the slit, since for any choice of a z , the 
factor 1131 immediately makes p{z\) 

zero outside the slit. However, inside the slit both the PY 
and the HLR equations take p Zl — p B and thus expand 
about the uniform bulk reference state, which gives a 
very poor description of the dilute 2D gas in the slit. 
As L s they predict a much lower limiting density2£ 
than given by eq[2S] These problems arise only at small 
separations of the order or smaller that of the hard sphere 
diameter d. At larger separations and in the one wall 
limit both theories give much more satisfactory results. 



VI. A NEW CRITERION FOR CHOOSING THE 
REFERENCE DENSITY IN THE SLR EQUATION 

The above discussion has shown the versatility of the 
SLR equation and its ability to give exact limiting results 
in several specific cases with proper choice of the p r . It 
has also shown that inaccuracies arise in some cases from 
the prescribed local choices made by the PY and HLR 
equations. Thus we need a more general and systematic 
way to choose p r in the SLR equation. To that end we 
first look more closely at the reasons why HLR fails in 
some cases. 



A. Limitations of the HLR prescription 

The HLR equation expands about the hydrostatic den- 
sity p T = p(p B — (fi(r)) . This depends on the external field 
too locally in cases where the external fields varies signifi- 
cantly in local regions much smaller than the correlation 
length of the fluid. In such cases, the proper density 
p r to expand about is often not the hydrostatic density, 
but a nonlocal extrapolation using the density of the sur- 
roundings, as illustrated by the tiny field and narrow slit 
examples discussed above. To use the SLR equation 1131 
to improve on the HLR approximation, we need a new 
way to choose p r that can account for this extrapolation 
of the local uniform system in such cases, while not spoil- 
ing the good results of the HLR equation in most other 
limits. We describe below our first attempt to develop 
such a computationally useful criterion. 
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B. The insensitivity criterion 

If all terms are exactly retained in eq 1121 it should 
be invariant with respect to a simultaneous shift of the 
chemical potential pP and the external field tfi. Thus 
cq ^] should hold for all choices of p r . However only 
certain choices of p r can efficiently truncate the series 
at low orders. One possible criterion for a truncation 
is to choose p r that minimize the contribution from 
the quadratic term Q(ri) in eq 1121 However, unlike 
c^(ri2',p), c( 3 )(ri,r2,r3;p) is often not available ana- 
lytically (and accurately) and the 6-dimensional integral 
of Q(r\) in eq. 0]is very computational demanding. 

To circumvent the difficulty of dealing with Q(r) di- 
rectly, a reasonable alternative is to consider how the 
predictions of the SLR equation change as p r is varied. 
Since the SLR equation is a truncation of the exact se- 
ries in eq 1121 it is certainly not invariant with respect to 
variation of any p r . However if the truncation is accurate 
for some particular choice of p r , then in effect the higher 
order terms in the series have been taken into account. 
Thus the SLR equation should be relatively insensitive 
to small variations about the particular p r that make the 
higher order corrections to the SLR equation small. This 
condition need not be exact, even in the case of a tiny 
field where the series truncates exactly, but it seems likely 
that it could produce reasonable choices for p r in many 
cases. 

This leads to the following self-consistent condition for 
the density given by the SLR equation [ED 

Sp( ri )/Sp" =0, Vn,r 2 , (27) 

expressing the insensitivity of the density with respect 
to variations in p r . Differentiating both sides of the 
SLR equation 1131 and collecting the terms arising from 
5p(ri)/Sp r2 (for details of the derivation, see the ap- 
pendix), the insensitivity criterion (IC) arising from eq 
E| can be written as: 

~ ri = /dr 2 Ty(|n-r 2 |;^)p(r 2 ) 

9 fdv 2 W(\r 1 -r 2 \;^) " 1 SJ 

where 

W{r 12 -p) =c«(p)c (2) (n 2 ;p) + c (2) (r 12 ;p) (29) 

and 

c«(p) =dc^(p)/dp; c {2 \r l2] p) = dc^{r 12 -p)/dp. 

(30) 

Because the function VK( |ri — r 2 |; p) in eqEHIhas range of 
c( 2 )(|ri — r 2 |; p), ea 1281 shows that p ri can be interpreted 
as the full density p(r) averaged over the range of the 
fluid's correlation length around the point ri , using a self- 
consistent weighting function W that itself depends on 
p ri . We will refer to the resulting p r as the smoothed ref- 
erence density in what follows. Some versions of weighted 
density functional theory have used similar weighted den- 
sities, though the detailed implementation and justifica- 
tion are rather different. 4,14 



Equation [35] derived from the IC can then be solved 
along with the SLR equation to determine both the full 
density p(r) and the smoothed density p r . We refer to 
these coupled equations as the IC equations. The IC 
equations can be solved numerically by iteration with the 
same methods used to solve the PY or HLR equations. 



C. Behavior of IC equations in limiting cases 

We first verify that the IC equations can give accurate 
results in limiting cases where the proper choice of p T is 
known. In the hydrostatic limit where the external field 
is very slowly varying, p(r) will reduce to the hydrostatic 
density p r = p(p B — ^(r)), as given by the HLR equation. 
In this same limit p(r 2 ) in the IC eauationl28lcan be ap- 
proximated by p(r±) and taken outside the integral. This 
gives p r = p(r) and hence p r = p r . The IC equations thus 
reduce to the HLR equation for slowly varying fields and 
recover the hydrostatic limit correctly. However, because 
of the averaging in eq|23 in other limits the IC choice of 
p r is less local than the HLR choice p r and tends to smear 
out the nonuniformity caused by external perturbations 
in small regions. 

For the narrow slit limit discussed above, the p r given 
by eq |2HI correctly approaches zero as L — > 0, since the 
Boltzmann factor e _/3< ^ r ' ensures that p(r) is zero inside 
the walls. Thus the IC equations are exact in the narrow 
slit limit as L s — > and correct the poor predictions of 
both the HLR and PY equations. 

For tiny fields, the IC choice in eq |2HI strictly repro- 
duces the exact bulk density choice p r = p(p B — 4>(oo)) 
only in the limit where V* — ► 0, and is not exact for 
general tiny fields. However p v is generally very close to 
that of the bulk density because the tiny field volume 
has little weight in the averaging. Thus the IC equations 
can be expected to agree reasonably well with exact re- 
sults for general tiny fields, as will be shown in the next 
section. 



VII. NUMERICAL RESULTS 

We first consider the density response of a hard sphere 
fluid at a moderate packing fraction r\ = 0.3 to a series 
of spherical model potentials of varying ranges and dif- 
ferent signs. In particular we studied repulsive (attrac- 
tive) step functions of height 3/csT (~3ksT) and dif- 
ferent ranges and "triangular" fields that start with the 
same height at the origin and vary linearly in r to the cut- 
off. Hard sphere cavities with the same cutoffs were also 
studied. Numerical solutions of the IC, HLR, and PY 
equations are presented, together with results of Monte 
Carlo simulations for the fluid's density response. The 
PY bulk direct correlation function (ri 2 ; p), which is 
very accurate at these densities, was used in the theoret- 
ical calculations^ The Carnahan-Starling 21 equation of 
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Tiny spherical repulsive potentials, r|=0.3 
Width=0.3d Width=0.5d 



Tiny spherical attractive potentials, t|=0.3 
Width=0.3d Width=0.5d 
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FIG. 1: The density response to tiny repulsive external fields 
of different widths W is plotted. The curves are the pre- 
dictions by the IC, PY and HLR approximations. The PY 
approximation is exact for the density values inside the non- 
zero field region. All the external fields are spherical. "Step" 
denotes a spherical step function where f3cf)(r) — 3, Vr < W 
and P4>(r) = 0, otherwise. "Triangle" refers to the poten- 
tial j3(j>(r) — 3 — 3r/W that has the same height as the 
step potentials but decays linearly to zero at r = W with 
/3(f>(r) = 0, otherwise. "Cavity" refers to the hard core po- 
tential P4>{r) — oo,Vr < W and P4>(r) = 0, otherwise. The 
form of the potentials /3<f>(r) are illustrated in the insets. For 
the cavity potentials, the PY and HLR approximations give 
identical density solutions. The bulk fluid's packing fraction 
r] = Tvp B d 3 /6 is r] = 0.3. 



state was used for computing the density p r = p(fi B —a r ) 
of the locally shifted uniform system. 



A. Tiny fields 

For tiny fields, all results should be compared to the 
PY approximation, which is exact for such fields (sub- 
ject only to the very small errors in the PY bulk direct 
correlation function). As can be seen in Figures ^ and 
13 the HLR equation is exact only for tiny cavities. For 
finite tiny fields, its major errors occur in the tiny region 
where the field is nonvanishing and rapidly varying. The 
density response predicted by the HLR equation often 
exhibits a negative region where the external field varies 
most rapidly. The IC approximation, on the other hand, 
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FIG. 2: The curves plotted here follow the same legend con- 
ventions used in Figure and are computed at the same bulk 
packing fraction r\ — 0.3. The external fields calculated here 
are all tiny spherical attractive potentials. Step and triangle 
potentials are defined as in Figure except that the sign of 
the potential is negative. 



Extended spherical repulsive potentials, T|=0.3 
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FIG. 3: Density response to extended spherical repulsive po- 
tentials of varying widths. Conventions and bulk packing frac- 
tion are those of Figures Q and |5] 



n 



Extended spherical attractive potentials, r|=0.3 
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FIG. 4: Density response to extended spherical attractive po- 
tentials of varying widths. Conventions and bulk packing frac- 
tion are those of Figures and [5] 



in general agrees with the PY approximation much bet- 
ter in the tiny field region and in particular eliminates 
the negative densities given by HLR. However the IC 
equations are not exact for tiny fields, and tend to over- 
estimate the contact densities. 



B. Extended fields 

For more extended fields, none of the approximations 
are exact, so Monte Carlo simulations were carried out 
to test the various theories. As shown in Figures 02 and 
QJ as the range of the step and the triangle potentials 
becomes wider, the HLR approximation becomes more 
accurate. However, it can still exhibit negative densities 
in its solution for both repulsive and attractive step po- 
tentials, especially for the narrower steps. For the same 
width of the potentials, the HLR equation does better 
in predicting the response to the triangle potentials than 
to the step potentials. This agrees with the expectation 
that the HLR equation should be more accurate when 
the external field is more slowly varying. 

The PY approximation, on the contrary, becomes less 
accurate when the field's width increases, as is seen in 
Figure |3J This is because the PY expansion about the 
bulk density and extrapolation into the region where the 
external field is non-vanishing becomes less and less jus- 
tified when the range of the potential increases. This 
problem with the PY approximation becomes much more 
acute for attractive potentials, where its errors are mag- 
nified by the large Boltzmann factor in eq 03 and the 
results are so poor that we do not show them in Figure 
0] Indeed, the PY approximation for nonuniform fluid 
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FIG. 5: The left graph gives the bowl potentials /3(f)(r) = 
A{r 2 - 4), with A = 1 for the shallow bowl, A = 2 for the 
deeper bowl. The right graph compares the hydrostatic den- 
sities p T of the external fields to the p T given by the IC equa- 
tions. 



is hardly ever applied in practice except for strongly re- 
pulsive potentials, where the value of p r in the repulsive 
region is essentially irrelevant. 

The IC approximation is able again to correct the nega- 
tive densities given by HLR and, most notably, to capture 
the highly nontrivial density profile inside both positive 
and negative step potentials due to the packing of the 
hard spheres. However, for spherical cavities, although 
it is known that the PY (and HLR) approximation con- 
sistently predicts a contact density lower than the ex- 
act valuet^iii the IC noticeable over-corrects the contact 
densities. This has a deleterious effect on the rest of the 
profile, especially near the cavity region. 

In Figures|S]and|n|we consider the density response of a 
hard sphere fluid with bulk packing fraction rj = 0.25656 
to two deep attractive parabolic bowl potentials. Figure |3 
shows the bowl potentials on the left and the reference 
densities p T and p r for the HLR and IC equations on the 
right. Results for the HLR and IC approximations are 
compared to MC simulations 10 in Figure For the shal- 
lower bowl potential, both the IC and HLR approxima- 
tions agree well with the Monte Carlo simulations. The 
reference density p r of the IC approximation is very close 
to that of the hydrostatic density p T , as can be seen in the 
right graph of Figure except that p r varies smoothly 
near the edge of the bowl, while the hydrostatic density 
has a discontinuous derivative. 

For the deeper bowl potential, both approxima- 
tions deviate noticeably from the simulation data, but 
nonetheless capture the nontrivial oscillatory density pro- 
file inside the bowl. In particular, both reproduce the 
density minimum at the center of the bowl, where the 
external field is actually most attractive, due to nonlocal 
effects from packing of the hard spheres. However, the 
HLR density becomes negative at the bottom of the bowl, 
while the IC density remains positive, though somewhat 
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FIG. 6: Density responses to attractive parabolic bowl poten- 
tials. Here n = 0.25656. 



A=10, o/d=0.5 



A=10, o/d=1.5 



P(r) 











! — |C 






■■■■ HLR 






o MC 




I 


— WW 






11=0.3 










r/d 



FIG. 7: The density response to two external fields of the 
form Aerfc(r / a) / r . Note that the left vertical axis refers to 
the density curves while the right vertical axis labels f3<f)(r). 
Here r\ — 0.3. 



lower than the MC result. The reference density p r for 
the IC method for the deeper bowl potential has more 
oscillations than that of the shallow bowl potential, and 
exhibits a maximum at the center, which is the key for 
keeping its predicted full density positive. Once outside 
the bowl, all approximations agree well with the MC re- 
sult. 

Figure plots the density response to a soft continu- 
ous repulsive potential of the form (3(j>(r) = Aeiic(r / a) / r. 
This potential is important in our theory of ionic fluids, 22 
but for our purposes here just serves as an example of a 
soft repulsive potential. Here all the approximate results 
agree quite well with the simulations, except that the 
HLR equation again shows a narrow negative density re- 
gion for the more rapidly varying potential (left graph of 
Figure EJ . 

Finally, in Figure the density response to a planar 
hard wall determined using the IC and HLR approxima- 
tions is compared to the results of the generalized mean 
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FIG. 8: Density response to one planar hard wall positioned 
at z = 0. The bulk packing fraction is r\ — 0.314. The "ref" 
curves give the p z used with the different approximations. 
The p z for the GMSA approximation is obtained by using the 
GMSA p(z) in the IC equation and solving inversely for the 
corresponding p z . 



spherical approximation (GMSA). 23 The GMSA fits the 
contact density at the wall to the bulk equation of state 
using an exact sum rule, and is known to be very accu- 
rate for such systems. Thus it can be used as a bench- 
mark for the other approximations. As is well known, 
the HLR approximation (equivalent in this case to the 
PY approximation) agrees well with the GMSA except 
for its consistent underestimate of the contact density. 

The IC result is now significantly worse than it was in 
the other cases. The main problem arises from a severe 
overestimate of the contact density, which goes on to spoil 
the rest of the density profile. The state shown is at a 
moderate packing fraction r\ = 0.314 and the errors get 
even larger at higher densities. This problem is similar 
to that seen for the cavity potentials presented before. 

The IC method can also be used in an inverse way to 
determine what p r is needed to obtain a given density 
p(r) as a solution. We determined the p r associated with 
the accurate GMSA p(r) in this way, and it can be seen 
in Figure |S| that this GMSA p r deviates from the bulk 
density (used in the HLR equation) only very near the 
wall. The IC p r shows a similar deviation, but has more 
oscillations away from the wall and dips less low near the 
wall. 

This example shows that relatively small changes in 
p r near the wall can have large effects on the predicted 
density profile near the wall in the SLR equation. The 
fact that the HLR choice, clearly appropriate for slowly 
varying fields, continues to give reasonably good results 
for single hard walls and hard cavities seems somewhat 
fortuitous, as illustrated by the errors HLR produces for 
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rapidly varying but finite repulsive potentials. See, e.g., 
the left graph in Figure [7| The IC method, which gives 
very good results in most other limits, correctly predicts 
a positive correction to the HLR/PY contact density but 
overestimates its magnitude. Only a small change in p r 
would be needed to produce very good results, but the 
IC prescription is not able to determine this in advance. 

VIII. FINAL REMARKS 

The SLR equation provides a versatile framework for 
computing the density response of hard sphere fluids to 
general external fields. Since it satisfies the exact shift- 
ing property of the grand ensemble, it can accurately 
describe two important limits: slowly varying external 
fields (hydrostatic limit) and the opposite limit where 
field can be very rapidly varying but only perturbs the 
fluid in tiny or narrow regions. Errors in specific pre- 
dictions arising from the linear truncation in the SLR 
equation in other cases can be minimized by expanding 
about a different uniform system at every point in space. 

In practice there could be different prescriptions for 
how these local uniform systems could best be chosen 
in particular applications and for specific properties. In 
many cases the simple local HLR choice is quite sufficient. 
However this has problems for rapidly varying but finite 
fields and the SLR equation allows other choices. A gen- 
eral idea often used in other expansions in liquid state 
theory is to choose a reference density p r that at least 
minimizes the quadratic correction to the SLR equation 
in eq However, this is very complicated, and there 
will still be unknown contributions from the higher order 
terms. 

Here, as a first attempt, we have devised a compu- 
tationally efficient insensitivity criterion (IC), based on 
the idea that the SLR equation should be insensitive to 
small variations of p r . This property would be exactly 
satisfied if all terms in the expansion were taken into ac- 
count, and by imposing it self consistently on the SLR 
equation we hope to generate a truncation where the con- 
tributions from higher order terms are indeed small. The 
resulting IC method is extremely successful in correcting 
the negative density regions that the HLR often exhibits 
for rapidly varying finite fields, and although not exact, 
the IC method also shows considerable improvement over 
HLR for the tiny field cases. Moreover, the IC method is 



exact for narrow slits as the slit width L s — > 0, while the 
HLR and PY approximations have very significant errors 
in this limit. 

However, the IC method tends to overestimate the con- 
tact value of the density response to single hard core cav- 
ities of all sizes, and this damages the accuracy of the rest 
of the density profile. In practice this need not be a sig- 
nificant limitation, since these cases are reasonably well 
treated by the HLR and PY approximations. Other spe- 
cific conditions for hard core potentials, such as the sum 
rule used in the GMSA approximation, could be taken 
into account to improve the IC method in this limit. 

But it seems conceptually worthwhile to try to choose 
p r more generally within the SLR framework so that ac- 
curate results naturally arise in this limit as well. To that 
end we believe it would be profitable to further study the 
tiny field limit, where similar problems are encountered, 
to gain additional insights into the optimal choice of the 
reference density for the SLR equation. We also need 
more information about the analytic nature and unique- 
ness of solutions the IC or similar criteria can provide. 
We have solved the resulting nonlinear equations numer- 
ically by iteration and have found a stable self consistent 
solution. But there could be other solutions, or other 
branches that only a small change in the IC could favor. 
We have preliminary evidence that in the hard wall limit 
an alternative branch may exist for high density states 
that could give much better results, and plan further re- 
search along these lines. 
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APPENDIX: INSENSITIVITY CONDITION FOR 
THE SLR EQUATION 

After carrying out the functional derivative 
Sp(ri)/5p r2 in eq \TJ\ on the p{ri) given by the SLR 
equation E| one can rewrite the resulting expression as 



dr 3 



S( Tl -T 3 )-p^e-^ ri ^c^(r 13 ;p^) 



M^) 
dp** 



p^ e -^ r Hri)s(r 1 -r 2 )\c^(p r ^)- / dr 3 c( 2 )(r 13 ; p ri ) 



+ / dr 3 \c^(p^)c^(r 13 ;p^) + c^(r 13 -,p^)\ (p(r 3 )-p ri ) 



(A.l) 



12 



A special case of eq [S] relating the n — 1 and nth order direct correlation functions can be written as 

c (n_1) (ri, ...,r„_i;p) = / dr„c (n) (ri, r n _i, r„; p). 



(A.2) 



Using this on the right hand side of ea IA.1I the first two terms in the curly brackets cancel. Thus, the requirement 
that <5p(r 3 )/<5/5 r2 in eq. I A. II vanish for all r 2 and r3 as in eq 123 then implies that 



dr 3 



c(D ) c (2) ( ri3 ; ^ ) + c(2) (ri3 . ^ ) (p( ra ) - ^ ) = 0, 



(A.3) 



for all ri, from which ea 1281 follows. Note that by using eq IA.2I the equation above can be written as 

dv 2 J dr 3 [c( 2 )(r 12 ;p ri )c^(r 1 3;p ri ) + c( 3 )(r 1 ,r 2 ,r 3 ;p ri )] (p(r 8 ) - p ri ) = (A.4) 
This could also be derived by making the following approximation for Q(r) in eqQ] 

QfrOwMn)-^) y dr 2 y dr 3 [c (2) (^i 2 ; /5 ri )c< 2 > (r 13 ; p ri ) + c (3) (n, r 2 , r 3 ; ^ )] (p(r 3 )-p ri ), (A.5) 

I 

i.e., by assuming that p(r 2 ) differs little from p(ri) in Requiring that this approximation for the quadratic term 
the region of integration near ri in the definition for Q. vanish then gives ea I A. 41 
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